//Nathan Popham
//4003-573-01: Procedural Shading
//Lab 2
//
//Phong shader
//
//


surface
realPhong ( float Ka = 0.2, 
            Kd = 0.5, 
            Ks = 0.5,
            roughness = 0.1; 
	        color specularcolor = (1, 1, 1);
          )
{

    float sqr (float x) { return x * x; }


    color
    LocIllumWardAnisotropic (normal N;  vector V;
                             vector xdir;  float xroughness, yroughness;)
    {
        float sqr (float x) { return x*x; }

        float cos_theta_r = clamp (N.V, 0.0001, 1);
        vector X = xdir / xroughness;
        vector Y = (N ^ xdir) / yroughness;

        color C = 0;
        extern point P;
        illuminance (P, N, PI/2) {
        /* Must declare extern L & Cl because we're in a function */
        extern vector L;  extern color Cl; 
        float nonspec = 0;
        lightsource ("__nonspecular", nonspec);
        if (nonspec < 1) {
            vector LN = normalize (L);
            float cos_theta_i = LN . N;
            if (cos_theta_i > 0.0) {
            vector H = normalize (V + LN);
            float rho = exp (-2 * (sqr(X.H) + sqr(Y.H)) / (1 + H.N))
                / sqrt (cos_theta_i * cos_theta_r);
            C += Cl * ((1-nonspec) * cos_theta_i * rho);
            }
        }
        }
        return C / (4 * xroughness * yroughness);
    }

    
    float diff = 0;
    float spec = 0;
    normal Nn = normalize(N); 
    
    extern point P;
    illuminance (P, Nn, PI) {
       extern vector L;
       normal Ln = normalize(L);
       float NdotL = max(Nn . Ln, 0.0);
       
       if(NdotL > 0.0)
       {
           diff += Kd*NdotL;
           normal R = normalize(reflect(-Ln, Nn));
           float RdotV = max(R . -normalize(I), 0.0); 
           spec += Ks*(RdotV)*roughness;
       }
    }
    Ci = Cs * (Ka*ambient() + diff) + 
	 LocIllumWardAnisotropic(N, -normalize(I), (N[1], -N[0], 0), 2.0, 2.0);
    Oi = Os;  
    Ci *= Oi;
}
